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Abstract. We report on the first measurement of the three-point function with the position- 
dependent eorrelation funetion from the SDSS-III Baryon Oscillation Spectroscopic Survey 
(BOSS) Data Release 10 CMASS sample. This new observable measures the correlation be¬ 
tween two-point functions of galaxy pairs within different subvolumes, .^(r, r^), where is 
the location of a subvolume, and the corresponding mean overdensities, J(ri;,). This correla¬ 
tion, which we call the “integrated three-point function”, iCir) = (^(r,r^)J(ri)), measures 
a three-point function of two short- and one long-wavelength modes, and is generated by 
nonlinear gravitational evolution and possibly also by the physics of inflation. The iCir) 
measured from the BOSS data lies within the scatter of those from the mock galaxy catalogs 
in redshift space, yielding a ten-percent-level determination of the amplitude of i(^{r). The 
tree-level perturbation theory in redshift space predicts how this amplitude depends on the 
linear and quadratic nonlinear galaxy bias parameters ( 6 i and 62 ), as well as on the ampli¬ 
tude and linear growth rate of matter fluctuations {a^ and /). Combining i(^{r) with the 
constraints on bicr^ and fag from the global two-point correlation function and that on ag 
from the weak leasing signal of BOSS galaxies, we measure 62 = 0.41 ± 0.41 ( 68 % C.L.) 
assuming standard perturbation theory at the tree level and the local bias model. 
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1 Introduction 

Mode coupling plays a fundamental role in cosmology. A long-wavelength scalar density 
fluctuation modifies the formation of small-scale structure via gravitational evolution (see [1] 
for a review), and possibly also through the physics of inflation. This effect of long-wavelength 
modes manifests itself through a dependence of observables, e.g., the n-point statistics and 
the halo mass function, on the local long-wavelength overdensity, or equivalently, the position 
in space (see [2-4] for the n-point statistics and [5, 6] for the mass function). Measurements 
of spatially-varying, “position-dependent” observables capture the effects of mode coupling, 
and can be used to test our understanding of gravity and the physics of inflation. A similar 
idea of measuring the shift of the peak position of the baryonic acoustic oscillation in different 
environments has been studied in ref. [7]. 
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In this paper, we focus on the position-dependent two-point function. Consider a galaxy 
redshift survey. Instead of measuring the two-point function of galaxy pairs within the en¬ 
tire survey volume, we divide the survey volume into many subvolumes, within which we 
measure the two-point function of galaxy pairs. These two-point functions vary spatially 
from subvolume to subvolume, and the variation is correlated with the mean overdensities 
of the subvolumes with respect to the entire survey volume. As we show later in detail, this 
correlation measures an integral of the three-point function, which represents the response of 
the small-scale clustering of galaxies (as measured by the position-dependent two-point func¬ 
tions) to the long-wavelength density perturbation (as measured by the mean overdensities 
of the subvolumes) [2]. 

Not only is the position-dependent correlation function conceptually straightforward 
to interpret, but the computational requirement for measuring three-point statistics is also 
largely alleviated. The usual three-point correlation function measurements rely on finding 
particle triplets with the naive algorithm scaling as where Apar is the number of par¬ 
ticles. Current galaxy redshift surveys contain roughly a million galaxies, and we need 50 
times as many random samples as the galaxies for characterizing the survey window func¬ 
tion accurately (see, e.g., [8]). Counting triplets thus becomes computationally challenging. 
Similarly, the measurement of the three-point function in Fourier space, i.e., the bispectrum, 
requires counting of all possible triangle configurations formed by different Fourier modes, 
which is also computationally demanding. This explains why only few measurements of the 
three-point function of the large-scale structure have been reported in the literature [9-17]. 

The computational requirement is alleviated for the position-dependent correlation func¬ 
tion technique because we explore a subset of triplets corresponding to the “squeezed config¬ 
urations” of the three-point function; namely, two short-wavelength modes correlated with 
one long-wavelength mode. We only need to count particle pairs for measuring the two-point 
function in subvolumes, which scales as The scaling is further improved because the 

number of particles in each subvolume decreases by the number of the subvolumes A*; hence, 
it scales as As(Apar/As)^ = Np^^/Ng. The position-dependent correlation function technique 
is thus particularly efficient for extracting the information of the squeezed-limit bispectrum, 
which we shall demonstrate in appendix D, while it is relatively insensitive to the bispectrum 
in other configurations. 

In this paper, we report on the first measurement of the three-point function with the 
position-dependent correlation function from the SDSS-III Baryon Oscillation Spectroscopic 
Survey Data Release 10 (hereafter BOSS DRIO) CMASS sample [18, 19]. We compare this 
measurement with those from the PTHalos mock catalogs [20-22]. While the mocks were 
designed to reproduce the global two-point function of the BOSS DRIO CMASS sample, 
it is not guaranteed that they can reproduce the three-point function as measured by the 
position-dependent correlation function. We shall show that the position-dependent corre¬ 
lation functions from the real data and the mocks are consistent with each other. Finally, 
we use tree-level perturbation theory to predict the position-dependent correlation function 
as a function of the galaxy bias parameters and the cosmological parameters, and determine 
the quadratic nonlinear bias parameter of the BOSS DRIO CMASS sample by combining the 
constraints from the position-dependent correlation function, the global two-point function, 
and the weak lensing signal. 

The rest of the paper is organized as follows. In section 2, we define the position- 
dependent correlation function and the integrated three-point function, and describe the 
tree-level perturbation theory prediction for the integrated three-point function in redshift 
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space. In section 3 and 4, we apply the position-dependent correlation function technique to 
the mocks and the BOSS DRIO CMASS sample, respectively. The cosmological interpretation 
of the measurements is given in section 5. We conclude in section 6. In appendix A, we test 
our estimator using Gaussian realizations. In appendix B, we study the effects of using 
extended models of the bispectrum with the effective F 2 and G 2 kernels and a tidal bias. In 
appendix C, we compare the mocks and BOSS DRIO CMASS samples in different redshift 
bins. In appendix D, we use the Fisher matrix calculation to demonstrate the information 
content of the position-dependent correlation function. Throughout the paper we adopt 
the cosmology of the mocks as our fiducial cosmology, i.e., a flat ACDM cosmology with 
Dm = 0.274, = 0.0224, h = 0.7, ag = 0.8, and Ug = 0.95. 

2 Position-dependent correlation function and the integrated three-point 
function in redshift space 

2.1 Position-dependent correlation fnnction 

Consider a density fluctuation field, (5(r), in a survey (or simulation) volume W- The mean 
overdensity of this volume vanishes by construction, i.e., ^ ^ fy <5(r) = 0. The global 

two-point function is defined as 


^(r) = ((5(x)(i(x-Fr)) , 


( 2 . 1 ) 


where we assume that (5(r) is statistically homogeneous and isotropic, so ^(r) depends only 
on the separation r. As the ensemble average cannot be measured directly, we estimate the 
global two-point function as 

^ ^ (5(r + x)(5(x). (2.2) 

7’ t/ t/x.x~|“r^Vr 


The ensemble average of eq. (2.2) is not equal to ^(r). Specifically, 


I r d‘^r r 

VrJ 


m) 


x,x-|-rSVr 


d^x (<5(r + x)<5(x)) =^(r):^ / ^ 


d-^x . 


' XjX-TrgVr 


(2.3) 


The second integral in eq. (2.3) is I 4 only if r = 0, and the fact that it departs from I 4 is 
due to the finite boundary of R.. We shall quantify this boundary effect later in eq. (2.8). 

We now identify a subvolume Vy centered at r^, and compute the mean overdensity 
and the correlation function within Vy- The mean overdensity is 


<5(rL) = ^ [ d^r (5(r) 
Vl Jvy 


I 


d^r (5(r)IF(r — tl) , 


(2.4) 


where W (r) is the window function. Throughout this paper, we use a cubic window function 
given by 


3 

W{r) = WUr) = n 9{n), Bin) 

i=l 


1, \ri\ < T/2, 

0, otherwise , 


(2.5) 


where L is the side length of Vy. The results are not sensitive to the exact choice of the 
window function, provided that the separation between galaxy pairs is much smaller than L. 
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While (5 = 0, 5(ri) is non-zero in general. In other wor(is, if S{ri) is positive (negative), then 
this subvolume is overdense (underdense) with respect to the mean density in Vr- 

Using the same window function, we define the position-dependent correlation function 
in the subvolume centered at as 

/ 

x,r+xeVL 

^ ^ J + 

This is essentially an estimator for a local two-point function. In this paper we shall consider 
only the angle-averaged position-dependent correlation function (i.e., the monopole) defined 
by 

i{r,rL) = J ^ l(r,rL) = J ^ J (fx (5(r-hx)(5(x)IUL(r-bx-rL)lTL(x-rL) . (2.7) 


^(r,rL) = ^ 


(5(r -|- x)(5(x) 

x)(5(x)lUL(r -F X - rL)lUL(x - r^) . (2.6) 


Similarly to that of the global two-point function, the ensemble average of eq. (2.7) is 
not equal to |(r). Specifically, 

{^{r,rL)) = :^ J ^ J (fx ((5(r -F x)5(x))lUL(r-F X - rL)lUL(x - rz,) 

= ^{r)^ j ^ j WLir + x)WL{x') = C{r)fhndvyir) , (2.8) 

where /bndry(^) is the boundary effect due to the finite size of the subvolume. While 
/bndry(?’) = 1 for r = 0, the boundary effect becomes larger for larger separations. The 
boundary effect can be computed by the five-dimensional integral in eq. (2.8). Alternatively, 
it can be evaluated by the ratio of the number of the random particle pairs of a given sep¬ 
aration in a finite volume to the expected random particle pairs in the shell with the same 
separation in an infinite volume. We have evaluated /bndry(?') in both ways, and the results 
are in an excellent agreement. 

As the usual two-point function estimators based on pair counting (such as Landy-Szalay 
estimator which will be discussed in section 3.2) or grid counting (which will be discussed 
in appendix A) do not contain the boundary effect, when we compare the measurements to 
the model which is calculated based on eq. (2.7), we shall divide the model by /bndry(?’) to 
correct for the boundary effect. 


2.2 Integrated three-point function 

The correlation between .^(r, r^) and 5(rz) is given by 


(|(r,rz)(5(rz)) = 


1 


cff 

dvr 

£7 




d^xi j d?X2 ((5(r-I- xi)(5(xi)(5(x2)) 

X WL{-r + xi - rz)lTz(xi - rz)IUz,(x2 - r^) 

/ C(r + xi-F rz,xi-h rz,X 2 rz) 

X lUz(r + xi)IUz(xi)IUz(x2) , 
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(2.9) 


where C(ri,r 2 ,r 3 ) = ( 5 (ri)h(r 2 )( 5 (r 3 )) is the three-point correlation function. Because of 
the assumption of homogeneity and isotropy, the three-point function depends only on the 
separations |rj — rj| for i / j, and so is independent of r^. Furthermore, as 

the right-hand-side of eq. (2.9) is an integral of the three-point function, we will refer to this 
quantity as the “integrated three-point function,” = (^(r, ri)h(r£,)). 

iC(r) can be computed if (C(ri,r 2 ,r 3 ) is known. For example, standard perturbation 
theory (SPT) with the local bias model at the tree level in real space gives 

K{r) = ^?<SPT(r) 6 ? 62 <fe 2(0 1 (2-10) 


where *Cspt and are given below. Here, hi and 62 are the linear and quadratic (nonlinear) 
bias parameters, respectively. Because of the high dimensionality of the integral, we use the 
Monte Carlo integration routine in the GNU Scientific Library to numerically evaluate iC,{f)- 
The first term, iCsPT, is given by [23, 24] 


CsPT(ri,r2,r3) 


y6(n2)6(?’23) + ^^12,2z[^[{rl2)(t)'l{r2:i) + Cz(’’23)</>z(n2)] 

4 f ^ (/>((ri2)(/>[(r23) _ 6(ri2)(/);(^23) _ ^l{r2‘iWi{ri2) 
7 1 ri2ri3 r23 ri2 


I ^ 

+ ^ 12,23 


^i{ri2) + 


3</>)(n2) 

ri2 


^l{r 23 ) + 


3(/>z(^23) ' 1 

^23 J J 


-|- 2 cyclic , 


( 2 . 11 ) 


where ri 2 = |ri — r 2 |, 1112,23 is the cosine between ri 2 and r 23 ,' refers to the spatial derivative, 
and 

/ d r rJJc 

k‘^Pl{k)smc{kr), ()),(r) = J ^ Pi{k)smc{kr), (2.12) 

with Pi{k) being the linear matter power spectrum, and sinc(x) = sin(x)/x. The subscript I 
denotes the quantities in the linear regime. 

The second term, i(^b 2 i is the nonlinear local bias three-point function. Since halos 
(galaxies) are biased tracers of the underlying matter density held, the local bias prescription 
yields the density held of the biased tracers as h/i(r) = 6i(5m(r) -|- + •••) where hi and 

62 are the linear and nonlinear biases, respectively, and 5ra{'^) is the matter density held [25]. 
The nonlinear bias three-point function is then 


Cb 2 (ri,r 2 ,r 3 ) = 6(?’i2)6(?’23) + 2 cyclic . (2.13) 


Figure 1 shows the scale-dependencies of iCsPT and i(^h 2 at z = 0 with Pi{k) computed 
by CLASS [26]. We normalize by (t\i, where 

= {6i{rLf) = —J^ Pi{k)\WU\^)\^ (2.14) 

is the variance of the linear density held in the subvolume Vl ■ The choice of this normalization 
will become clear in section 2.4 where we discuss iC, in the squeezed limit, i.e., r L. We 
hnd that the scale-dependencies of *Cspt(^) and iCz) 2 (’’) are similar especially on small scales. 
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Figure 1: Normalized iCsPT (solid) and (dashed) for L = 100 h ^ Mpc (red), 
200 h~^ Mpc (green), and 300 h~^ Mpc (blue) at z = 0. 


This is because the scale-dependence of the bispectrum in the squeezed limit is (see e.g. the 
appendix of ref. [2]) 


-Bspt —^ 


68 


1 din k^Pi{k) 
3 dink 


m)Pi{q ), 


Bb2 ‘^Pl{k)Pi{q) , 


(2.15) 


where k and q are the short- and long-wavelength modes, respectively. For a power-law power 
spectrum without features, the squeezed-limit i?spT and have exactly the same scale 
dependence and cannot be distinguished. This results in a significant residual degeneracy 
between bi and b 2 , and will be discussed in section 5. When r is small, i(^{r)/a\i becomes 
independent of the subvolume size. We derive this feature when we discuss the squeezed limit 
in section 2.4. 


2.3 Connection to the integrated bispectrum 

Fourier transforming the density fields, the integrated three-point function can be written as 


= 772 


d^qi 


dz’qfs 


yl 


(27r)®(5£)(qi + q 2 + + q 2 + q4 + q5)<5D(q3 + qe) 


(27r)^ J (27r)3 

xil(qi,q2,q3)VFL(q4)VFi(q5)l^L(q6)e*'"('^^+'l^)-^^-(‘^^+‘^®+'^'^)] 


d?k 

(27r)3 


iB{k)e 


2rk 


(2.16) 
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where i?(qi,q 2 ,q 3 ) is the bispectrum of the tracers, and 


= j J S(k-qi,-k + qi+q3,-q3)VbL(qi)bbL(-qi-q3)bbL(q3), 

(2.17) 

is the integrated bispectrum as defined in eq. (2.7) of ref. [2]. Eq. (2.16) shows that the inte¬ 
grated three-point function is the Fourier transform of the integrated bispectrum. Similarly, 
the angle-averaged integrated three-point function is related to the angle-averaged integrated 
bispectrum, iB{k) = (47r)“^ f (Pk iB(k), as 

/ k'^rlk 

^ ^ iB{k) sinc(A:r). (2-18) 


2.4 Squeezed limit 

In the squeezed limit, where the separation of the position-dependent correlation function is 
much smaller than the size of the subvolume (r <C E), the integrated three-point function has 
a straightforward physical interpretation [2]. In this case, the mean density in the subvolume 
acts effectively as a constant “background” density. Consider the position-dependent corre¬ 
lation function, ^(r, r/,), measured in a subvolume with overdensity (5(r£,). If the overdensity 
is small, we may Taylor expand ^(r, r^) in orders of 5 as 




?(r)l5=o + 


d6 


6 + 0{P) . 

( 5=0 


(2.19) 


The integrated three-point function in the squeezed limit is then, at leading order in the 
variance {5'^) (dropping 5 = 0 in the subscript of the derivative term for clarity), given by 


iC(r) = («r,ri)«(ri)) = + 0(P) . (2,20) 

As {P) = fC(r) normalized by (t£ is d^{r)/d6 at leading order, which is the linear 
response of the correlation function to the overdensity. Note that in eq. (2.20) there is no 
dependence on the subvolume size apart from as shown also by the asymptotic behavior 
of the solid lines in figure 1 for r —)• 0. 

As iCir) is the Fourier transform of iB(k), the response of the correlation function, 
d^(r)/d6, is also the Fourier transform of the response of the power spectrum, dP{k)/d5. For 
example, we can calculate the response of the linear matter correlation function, d^i{r)/d6, 
by Fourier transforming dPi{k)/dS = [68/21 — {l/3)d\nk^Pi{k)/dliik]Pi{k) [2]. In figure 2, 
we compare the normalized fCsPT(^') with d^i{r)/d5. Due to the large dynamic range of 
the correlation function, we divide all the predictions by ^{r). As expected, the smaller 
the subvolume size, the smaller the r for fCsPT(^) to be close to [l/^i{r)][d^i{r)/d6], i.e., 
reaching the squeezed limit. Specifically, for 100 h~^ Mpc, 200 h~^ Mpc, and 300 h~^ Mpc 
subvolumes, the squeezed limit is reached to 10% level at r r-u 10 h ^ Mpc, 18 /i ^ Mpc, and 
25 h~^ Mpc, respectively. 

^If 5 = (5; then ;. But S can in principle be nonlinear or the mean overdensity of the biased tracers, 

so here we denote the variance to be (t|. 
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Figure 2: The linear response function \1 /^i{r)][d^i{r)/dS\ (black solid) and the normalized 
zCspt(?’) for L = 100 h~^ Mpc (red dotted), 200 h~^ Mpc (green dashed), and 300 h~^ Mpc 
(blue dot-dashed). The light and dark bands correspond to ±5% and ±10% of the predictions, 
respectively. 


2.5 Bispectrum in redshift space 

To model iC{r) in redshift space, we need a model for the bispectrum in redshift space. SPT 
in redshift space at the tree level predicts the halo bispectrum with local bias as [27] 


B2_tree-ievei(ki,k2,k3) = 2 (^ 2 (ki, k 2 )Zi(ki)Zi(k 2 )Pi(/ci)Pi(/C 2 ) ± 2 cyclic], (2.21) 

with 


Zi{ki) 


■^2(ki, k2) 


= {h + ff4), 

= 6iP2(ki,k2) ± //r^G2(ki,k2) ± 


ffik 


/P 

ki 


{h + f^i^) + + fk'i) 


+ 


h 
2 ’ 
( 2 . 22 ) 


where F 2 and G 2 are the standard kernels of SPT [1], / = din D/din a is the logarithmic 
growth rate, ^ = k ■ r\os, Hi ^ ki ■ rios, and k = ki ± k 2 . 

The integrated three-point function is the Fourier transform of the integrated bispectrum. 
Thus, we can evaluate iC,{r) by using eq. (2.21) in eq. (2.17) and averaging over the angle of 
k as in eq. (2.18). This operation requires a nine-dimensional integral. On the other hand, if 
we have an expression for the three-point function in configuration space, such as eq. (2.11), 
we can use eq. (2.9) to evaluate iCir), which requires an eight-dimensional integral. We do 
not always have an analytical expression for the three-point function in configuration space; 
thus, we in general need to perform the nine-dimensional integral to obtain iC{r) from the 
bispectrum. 
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Nevertheless, to check the precision of numerical integration, we compare the results 
from the eight-dimensional integral in eq. (2.9) with eq. (2.11), and the nine-dimensional 
integral in eq. (2.17) and eq. (2.18) with il(ki, k 2 , ka) = 2(^2(ki, k 2 )Pi(A:i)P/(/i: 2 )-|-2 cyclic]. 
As the latter gives a noisy result, we apply a Savitzky-Golay filter (with window size 9 and 
polynomial order 4) six times. We find that, on the scales of interest (30 h~^ Mpc < r < 
78 h~^ Mpc, which we will justify in section 3.3), both results are in agreement to within 
2%. We repeat the same test for (eq. (2.13)), finding a similar result. As the current 
uncertainty on the measured integrated correlation function presented in this paper is of 
order 10%, we conclude that our numerical integration yields sufficiently precise results. 


2.6 Shot noise 

If the density field is traced by discrete particles, 6d{r), then the three-point function contains 
a shot noise contribution given by 


('5d(ri)(5rf(r2)(5rf(r3)) 




((I(ri)5(r2)<5(r3)) 

+2 cyclic 

L nirs) 


Spjri - r2)<5D(ri - r3) 
h(r2)h(r3) 

(2.23) 


where n(r) is the mean number density of the discrete particles. The shot noise can be safely 
neglected for the three-point function because it only contributes when ri = r 2 , ri = r 3 , or 
r 2 = r 3 . On the other hand, the shot noise of the integrated three-point function can be 
computed by inserting eq. (2.23) into eq. (2.9), which yields 


iCshotir) = ar)^ I ^ I 


+ 


_n{x + r + rL) n(x-|-ri)_ 


WL(x + r)WL(x) , (2.24) 


where we have assumed r ^ 0. If we further assume that the mean number density is 
constant, then the shot noise of the integrated three-point function can be simplified as 


<shot(?’) = ‘^C{r)^^fhndry{r) ■ (2.25) 

For the measurements of PTHalos mock catalogs and the BOSS DRIO CMASS sample, the 
shot noise is subdominant (less than 7% of the total signal on the scales of interest). 


3 Application to PTHalos mock catalogs 

Before we apply the position-dependent correlation function technique to the real data, we 
apply it to the 600 PTHalos mock galaxy catalogs of the BOSS DRIO CMASS sample in 
the North Galactic Cap (NGC). From now on, we refer to the real and mock BOSS DRIO 
CMASS samples as the “observations” and “mocks”, respectively. 

We use the redshift range of 0.43 < z < 0.7, and each realization of mocks contains 
roughly 400,000 galaxies. We convert the positions of galaxies in RA, DEC, and redshift to 
comoving distances using the cosmological parameters of the mocks. The mocks have the 
same observational conditions as the observations, and we correct the observational system- 
atics by weighting each galaxy differently. Specifically, we upweight a galaxy if its nearest 
neighbor has a redshift failure (tc^f) or a missing redshift due to a close pair (tCcp). We fur¬ 
ther apply weights to correct for the correlation between the number density of the observed 
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galaxies and stellar density (tfstar) and seeing {wsee)- We apply the same weights as done in 
the analyses of the BOSS collaboration, namely FKP weighting, rcpKP = [l + -Pu)h(z)comp]“^ 
[28], where Pu, = 20000 h~^ Mpc^, and n{z) and “comp” are the expected galaxy number 
density and the survey completeness, respectively, provided in the catalogs. Therefore, each 
galaxy is weighted by icboss = {wcp + w^f - l)rCstar?nseeW'FKP- 

In this section, we present measurements from mocks in real space in section 3.3 and 
redshift space in section 3.4. The application to the CMASS DRIO sample is the subject of 
section 4. 


3.1 Dividing the subvolumes 

We use SDSSPix^ to pixelize the DRIO survey area. In short, at the lowest resolution (res=l) 
SDSSPix divides the sphere equally into Ux = 36 longitudinal slices across the hemisphere 
(at equator each slice is 10 degrees wide), and each slice is divided into Uy = 13 pieces along 
constant latitudes with equal area. Thus, for res=l there are nxXUy = 468 pixels. In general 
the total number of pixels is n(, x n'y = (res Ux) x (res Uy) = (res)^ x 468, and in this paper 
we shall set res=1024. After the pixelization, the object (a galaxy or a random sample) 
has the pixel number {ix,iy)- 

We use two different subvolume sizes. To cut the irregular survey volume into subvol¬ 
umes with roughly the same size, we first divide the random samples at all redshifts into 10 
and 20 slices across longitudes with similar numbers of random samples; we then divide the 
random samples in each slice into 5 and 10 segments across latitudes with similar numbers 
of random samples. Figure 3 shows the two resolutions of our subvolumes before the redshift 
cuts. (Note that this resolution is different from the resolution of SDSSPix, which we always 
set to res=1024.) Each colored pattern extends over the redshift direction. Finally, we di¬ 
vide the two resolutions into three (zcut = 0.5108, 0.5717) and five (^Jcut = 0.48710, 0.52235, 
0.55825, 0.60435) redshift bins. 

As a result, there are 150 and 1000 subvolumes for the low and high resolution config¬ 
urations, respectively. The sizes of the subvolumes are approximately = 220 h~^ Mpc 
and 120 h~^ Mpc, respectively^. The fractional differences between the numbers of the ran¬ 
dom samples in subvolumes for the low and high resolutions are within and 33 %; 

respectively. Since the number of random samples represents the effective volume, all subvol¬ 
umes at a given resolution have similar effective volumes. We assign galaxies into subvolumes 
following the division of random samples. 


3.2 Estimators in the subvolumes 


In the subvolume, we measure the mean overdensity with respect to the entire NGC, hj, 
and the position-dependent correlation function, ^i{r). The mean overdensity is estimated by 
comparing the total weighted galaxies to the expected number density given by the random 
samples, i.e.. 


a Wr^i 


Wr,i Wr^tot 


(3.1) 


^SDSSPix: http: //dls .physics .ucdavis. edu/'scranton/SDSSPix 

^The shapes of the subvolumes are not exactly cubes. For example, for the high resolution, the ratios of 
square root of the area to the depth, ^JL^Ly |are roughly 0.78, 1.42, 1.51, 1.28, and 0.71, from the lowest 
to the highest redshift bins. The results are not sensitive to the exact shape of the subvolumes, as long as the 
separation of the position-dependent correlation function that we are interested in is sufficiently smaller than 
F 3 ;, Ty, and Lz- 
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Figure 3: Division of random samples into subvolumes with two resolutions in the RA-DEC 
plane. Each colored pattern extends over the redshift direction. 


where Wg^i and Wr^i are the total weights (wboss) of galaxies and random samples in the 
subvolume, respectively, and Ng is the number of subvolumes. 

We use the Landy-Szalay estimator [29] to estimate the position-dependent correlation 
function as 


^ DDi{r, fl) f Er Wr,i? - Er 
RRi{r, /i) [J2g Wg^i]^ - J2g wli J 


DRi(r, m) (Er - Er ^r,i) ^ ^ 

RRi{r,n) Y^gWg^iY^^Wr,i 


(3.2) 

where DDi{r, /i), DRi{r, ^), and RRi{r, /r) are the weighted numbers of galaxy-galaxy, galaxy- 
random, and random-random pairs within the subvolume, respectively, and /r is the 
cosine between the line-of-sight vector and the vector connecting galaxy pairs (ri — r 2 ). The 
summations such as Wr,i and Yg denote the sum over all the random samples and 
galaxies within the subvolume, respectively. The angular average correlation function is 
then = fo dfJ. |LS,i(D/^)- 

Eq. (3.2) estimates the correlation function assuming that the density fluctuation is 
measured relative to the local mean. However, the position-dependent correlation function 
defined in section 2 uses the density fluctuation relative to the global mean. These two 
fluctuations can be related by (5giobai = (1 + ^)<5iocai + d with S = hiocai/hgiobai ~ 1- Thus, the 
position-dependent correlation function, ^i(r), is related to the Landy-Szalay estimator as 

?*(^) = (1 + ^i)^^LS,i(^) + . (3-3) 


To compute the average quantities over all subvolumes, we weight by Wr^i in the corre¬ 
sponding subvolume. For example, for a given variable gi in the subvolume, the average 
over all subvolumes, g, is defined by 




9 = 


^r,tot 




(3.4) 


2=1 


Since the number of random samples in each subvolume represents the effective volume, the 
average quantities are effective-volume weighted. Eq. (3.4) assures that the mean of the 
individual subvolume over densities is zero, 


5 


1 


Wr^tot 


Ns 

di'Wr,i 

2=1 


1 




-'Vs 

E 


2=1 


1 

a 


--1 = 0 


a 


(3.5) 
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We also confirm that ^{r) from eq. (3.3) agrees with the two-point function of all galaxies in 
the entire survey, on scales smaller than the subvolume size, integrated three-point function 
in the subvolume of size L as 


1 


'^r,tot 


2 = 1 


^i{r)Si - 2 ^(r) 


(1 + a) Er 

« Er W.icomp^ . 


2 

r,i 


-1 


^ n^.icomp^ j 


Wr 


(3.6) 

where the second term in the parentheses is the shot noise contribution, and Ur^i and comp^ j 
are the expected galaxy number density and the survey completeness, respectively, of the 
random samples. Similarly, we estimate the shot-noise-corrected variance of the fluctuations 
in the subvolumes of size L as 


Ns 


<^L = 


1 “ 
— E 

Wr.tot “ 
’ 2 = 1 


(1 -|-a) Ylr'^r,i f ^ 

« Er ^r,iComp^ iW^ - hr.,iCOmp^ j 


W. 




(3.7) 


where the second term in the parentheses is the shot noise contribution. We find that the 
shot noise is subdominant (less than 10 %) in both iC{r) and cr^. 


3.3 Measurements in real space 

Figure 4 shows the measurements of the two-point function ^(r) from the entire survey (top 
left) and the normalized integrated three-point functions (bottom panels), iQ{r)/a\^ for the 
subvolumes of two sizes (220 h~^ Mpc in the bottom-left and 120 h~^ Mpc in the bottom- 
right panels). The gray lines show individual realizations, while the dashed lines show the 
mean. 

We now fit models of ^(r) and iQ{r)/cj\ to the measurements in 30 h~^ Mpc < r < 
78 h~^ Mpc. We choose this fitting range because there are less galaxy pairs at larger 
separations due to the subvolume size, and the nonlinear effect becomes too large for our 
SPT predictions to be applicable at smaller separations. For the two-point function, we take 
the Fourier transform of [30] 

Pg{k) = b\[Pi{k)e~^ -I-^mc-Pmc(^)] , (3.8) 

where bi is the linear bias, Pi{k) is the linear power spectrum, Tmc is the mode coupling 
constant, and 

PMc(fc)=2 I ^p,(g)p,(|k-q|)[T2(q,k-q)]E (3.9) 

Hence, ^g{r) = 6 f[ 6 ,< 7 „(?’) + ^mc^mc(^’)] with 

= / (03 V {Mc(r) = / (^ ■ (3.10) 

We use a fixed value of = 20.644 h~‘^ Mpc^. Varying it has only small effect on the other 
fitted parameters. For the integrated three-point function, we use the SPT calculation 

iCgjr) _ 6PCsPT(r) ^2<&2(r) 1 . . 

/bndry(r) ’ ^ ^ 
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Figure 4: (Top left) ^(r) of the mocks in real space. The gray lines show individual real¬ 
izations, while the dashed line shows the mean. The black solid line shows the best-fitting 
model. (Top right) x^-histogram of the 600 mocks jointly htting the models to ^(r) and 
iQ{r)l(j\ in real space. The dashed line shows the y^-distribution with d.o.f.=36. (Bottom 
left) iC,{r')la\ of the mocks in real space for 220 h~^ Mpc subvolumes. (Bottom right) Same 
as the bottom left panel, but for 120 h~^ Mpc subvolumes. 


where iCsPT(^') and *Cfe 2 (^) computed from eq. (2.9) with eqs. (2.11) and (2.13), respec¬ 
tively, and is computed from eq. (2.14), using the subvolume sizes of L = 220 and 
120 h~^ Mpc and the redshift of z = 0.57. Note that the size of the subvolumes affects the 
values of cr| ^. We determine L by first measuring h\ using the real-space two-point function 
of the entire survey, and then find L such that assuming the cubic top-hat win¬ 

dow function'^. We find that these values {L = 220 and 120 h~^ Mpc) agree well with the 
cubic root of the total survey volume divided by the number of subvolumes, to within a few 
percent. 

We fit the models to ^(r) and iC,{r')la\ of both subvolumes simultaneously by minimizing 

= E ’ (3.12) 

ij 

where C~^ is the inverse covariance matrix computed from the 600 mocks, Di and Mj are the 
^In principle, the shape of the window function also affects cr|, but we ignore this small effect. 
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data and the model in the bin, respectively. The models contain three fitting parameters 
bi, 62 , and ^mc- 

The models computed with the mean of the best-fitting parameters of 600 mocks are 
shown as the black solid lines in figure 4. The best-fitting parameters are 61 = 1.971 ± 0.076, 
62 = 0.58 ± 0.31, and ^mc = 1-44 ± 0.93, where the error bars are 1-a standard deviations. 
The agreement between the models and the mocks is good, with a difference much smaller 
than the scatter among 600 mocks. Upon scrutinizing, the difference in ^(r) is larger for 
larger separations because the ht is dominated by the small separations with smaller error 
bars. On the other hand, for iC,(r^ja\ the agreement is good for both sizes of subvolumes 
at all scales of interest. This indicates that the SPT calculation is sufficient to capture the 
three-point function of the mocks in real space. 

The data points in hgure 4 are highly correlated. To quantify the quality of the fit, we 
compute the x^-histogram from 600 mocks, and compare it with the y^-distribution with the 
corresponding degrees of freedom (d.o.f.). There are 13 fitting points for each measurement 
(^(r) and two sizes of subvolumes for zC(r)/cr£) and three fitting parameters, so d.o.f.=36. 
The top right panel of hgure 4 shows the y^-histogram. The dashed line shows the x^- 
distribution with d.o.f.=36. The agreement is good, and we conclude that our models well 
describe both ^(r) and of the mocks in real space. 

Our 61 is in good agreement with the results presented in hgure 16 of ref. [16], whereas 
our 62 is smaller than theirs, which is ~ 0.95, by 1.2a. This may be due to the difference in the 
bispectrum models. While we restrict to the local bias model and the tree-level bispectrum, 
ref. [16] includes a non-local tidal bias [31-33] and uses more sophisticated bispectrum mod¬ 
eling using the effective F 2 kernel [34, 35]. In appendix B, we show that using the effective 
F 2 kernel and the non-local tidal bias in the model increases the value of b 2 , but the changes 
are well within the 1-a uncertainties. Also, the differences of the goodness of ht for various 
models are negligible. 

The htting range as well as the shapes of the bispectrum may also affect the results: 
the integrated correlation function is sensitive only to the squeezed conhgurations, whereas 
ref. [16] includes more equilateral and collapsed triangle conhgurations. Understanding this 
difference merits further investigations. 

3.4 Measurements in redshift space 

Figure 5 shows the measurements of ^(r) (top left) and i(^{r)/a\ (220 h~^ Mpc in the 
bottom-left and 120 h~^ Mpc in the bottom-right panels) of the mocks in redshift space. 
The gray lines show individual realizations, while the dashed lines show the mean. Similar 
to the analysis in section 3.3, we ht the models in redshift space to the measurements in 
30 h~^ Mpc < r < 78 h~^ Mpc. In this section, we use General Relativity to compute the 
growth rate, f{z) ^ which yields f{z = 0.57) = 0.751. We shall allow / to vary 

when interpreting the measurements in the actual data. 

Since there is no baryonic acoustic oscillation feature on the scales we are interested in, 
we model the redshift-space two-point correlation function as 

Cg,zir) = b\ [ 6 ,o-„(?’) + ^MC?Mc(?’)] Kp , (3.13) 

where ^i,(T„(r) and ^mc(^’) are given in eq. (3.10) and 

Rp = 1 + g/? + -^0^ 5 (3-14) 
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Figure 5: Same as figure 4 but in redshift space. 


is the Kaiser factor with /3 = f/bi [36]. As we do not include the subdominant term pro¬ 
portional to 62 in the two-point function, it only gives a constraint on bi, which we can then 
use to break the degeneracy with 62 in the integrated three-point function. We find that this 
simple modeling yields unbiased 61 and fulfills the demand. We calculate the redshift-space 
normalized integrated three-point function using SPT at the tree level, as described in sec¬ 
tion 2.5, and then correct for the boundary effect. The of the mocks in redshift space 
agrees with b\K ^ to percent level. The redshift-space models thus contain, as before in 
real space, the three fitting parameters, 61 , 62 , and Amc- We then simultaneously fit ^(r) and 
zC(r)/cr£ of both subvolumes by minimizing eq. (3.12). Figure 6 shows the correlation matrix 
{Cij in normalized by ^JCuCjj) estimated from the 600 mocks in redshift space. Because 
we normalize the integrated three-point function by the covariance between 
and is negligible. On the other hand, the covariances between and ^(r), between 

^(r), and between i(^{r')ja\ for two sizes of subvolumes are significant. 

The models computed with the mean of the best-fitting parameters of 600 mocks are 
shown as the thick solid lines in figure 5. The best-fitting parameters are bi = 1.931 ±0.077, 
62 = 0.54 ± 0.35, and Amc = 1-37 ± 0.82. The agreement between the models and the 
measurements in redshift space is as good as in real space. 

Again, our bi is in good agreement with the results presented in figure 16 of ref. [16], 
whereas our 62 is smaller than theirs, which is ~ 0.75, but still well within the 1-cr uncertainty. 
As noted in section 3.3, the adopted models of the bispectrum are different. In appendix B, 
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Figure 6 : Correlation matrix estimated from 600 mocks in redshift space. The figure shows 
a1 and iC,{r)/a\ of 220 h~^ Mpc subvolumes from bin 0 to 13, (j\ and iC,{r)l(j\ of 120 Mpc 
subvolumes from bin 14 to 27, and ^{r) from bin 28 to 40. 


we show that using the effective F 2 and G 2 kernels and the non-local tidal bias in the model 
increases the value of 62 - However, the changes are within the uncertainties, and the goodness 
of the ht is similar for different models. Thus, in this paper we shall primarily use the SPT 
at the tree level with local bias for simpler interpretation of the three-point function, but 
also report the results for the extended models. 

4 Measurements of the BOSS DRIO CMASS sample 

We now present measurements of the position-dependent correlation function from the BOSS 
DRIO CMASS sample^ in NGC. The detailed description of the observations can be found 
in refs. [18, 19]. Briefly, the sample contains 392,372 galaxies over 4,892 deg^ in the red- 
shift range of 0.43 < z < 0.7, which corresponds to the comoving volume of approximately 
2 h~^ Gpc^. We also weight the galaxies by tCBOSS to correct for the observational sys- 
tematics. We follow section 3.1 to divide the observations into subvolumes. However, the 
observations have their own set of random samples, which are different from the ones of the 
mocks (the random samples of the mocks have slightly higher n and different n(z)), so we 
adjust the redshift cuts to be Zcut = 0.5108, 0.5717 and Zcut = 0.48710, 0.52235, 0.55825, 
0.60435 for the two resolutions, respectively. The resulting properties of subvolumes of the 
observations and mocks are similar. 

The mocks are constructed to match the two-point function of the observed galaxies, but 
not for the three-point function. Hence there is no guarantee that the three-point function 
of mocks agrees with the observations. We can test this using our measurements. 

The measurements of ^(r) and iC,{r')la\ from the observations are shown as the solid 
lines in figure 7; the measurements of a\ is summarized in table 1. The measurements are 
consistent visually with the mocks within the scatter of the mocks®, and we shall quantify 

^Catalogs of galaxies and the random samples can be found in http://www.sdss3.org. 

®These measurements of i(,(r)la\ are done for one effective redshift. We compare of the ob¬ 

servations and mocks in different redshift bins in appendix C, finding that the observations and mocks are 
consistent at all redshift bins to within the scatter of the mocks. 
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Figure 7: Measurements of the BOSS DRIO CMASS sample (black solid lines). The gray 
lines show individual mocks in redshift space and the dashed line shows the mean of mocks. 
(Top left) ^(r), (Bottom left) iC,{r)l(y\ for 220 h~^ Mpc subvolumes, and (Bottom right) 
for 120 h~^ Mpc subvolumes. (Top right) y^-histogram of the 600 mocks jointly 
fitting the three amplitudes to ^(r) and iQ{r)/a\ in redshift space. The dashed line shows 
the x^-distribution with d.o.f.=38. The x^ value measured from the BOSS DRIO CMASS 
sample is 46.4. 




220 h-^ Mpc 
120 h~^ Mpc 

4.6 X 10“^ 5.6 X 10-'^ 4.9 x 10"^ 
2.4 X 10"2 1.3 X 10"3 2.5 X 10"2 


Table 1: Measurements of a\ of the mock catalogs and the BOSS DRIO CMASS sample. 


the goodness of fit using x^ statistics later. 

To quantify statistical significance of the detection of iC,{r)ja\ and the goodness of fit, 
we use the mean of the mocks as the model (instead of the model based on perturbation theory 
used in section 3.4), and fit only the amplitudes of iC,{r)la\, ^(r), and to the observations 
and the 600 mocks by minimizing eq. (3.12). Specifically, we use Oi{r) = Ai 0™°^’^(r) as the 
model, where Oi{r) = iQ{r)/a\, 02(r) = ^(r), and with the amplitudes Ai, A 2 , A 3 . 

Table 2 summarizes the fitted amplitudes. The l-cr uncertainties and the correlations 
are estimated from the 600 mocks. Since we normalize iC(^) by ct£, the correlation between 
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Ai A 2 A 3 

1 -a error 
best-fit (DRIO) 

0.12 0.03 0.04 

0.89 1.02 1.08 



(%l,^ 2 ) (%1,%3) (^2,%3) 

corr 

0.34 0.09 0.36 


Table 2: Results of fitting the amplitudes: Ai is iC,{r)jo\^ is ^(r), and ^3 is a\. (Left) 
The \-a uncertainties of the amplitudes estimated from the mocks, and the best-fitting 
amplitudes of BOSS DRIO CMASS sample with respect to the mean of the mocks. (Right) 
The correlation coefficients of the amplitudes. 


Ai and A 3 is small. On the other hand, A 2 and A 3 are correlated significantly because is 
an integral of the two-point function [eq. (2.14)]. 

Comparing the BOSS DRIO CMASS sample to the mean of the mocks, we find that 
zC(r)/cr£ is l-cr lower, ^(r) is unbiased (by construction of the mocks), and (t| is 2 -a higher. 
The result of Ai for the data is driven by the correlation between different separations of 
iC,(r)la\. On the other hand, the result of A 3 is driven by the positive correlation between 
^(r) and (j\. While of the data for two subvolumes are larger than that of the mocks 
but still at the boundary of the variances (see table 1 ), it requires an even higher A 3 to 
minimize when we jointly fit the three amplitudes. The fact that A 3 is larger than A 2 
is also possibly due the contributions to a\ from small separations (including stochasticity 
at zero separations), where the mocks were not optimized. We find Ai = 0.89 ± 0.12, i.e., a 
7Aa detection of the integrated three-point function of the BOSS DRIO CMASS sample. 

In order to assess the goodness of fit, we use the distribution of a histogram of which 
is shown in the top right panel of figure 7. In total there are 41 fitting points (13 fitting points 
for ^(r) and two sizes of subvolumes for iC^{r)ja\, and two fitting points for a\) with three 
fitting parameters, so d.o.f.=38. The x^ value of the observations is 46.4, and the probability 
to exceed this x^ value is more than 16%. Given the fact that the mocks are constructed to 
match only the two-point function of the observations, this level of agreement for both the 
two-point and integrated three-point correlation functions is satisfactory. 

5 Cosmological interpretation of the integrated three-point fnnction 

What can we learn from the measured In section 3.4, we show that the prediction 

for iQ (r) /based on SPT at the tree-level in redshift space provides an adequate fit to the 
mocks to within the scatter of the mocks; thus, we can use this prediction to infer cosmology 
from iQ{r)/a\. Note that any unmodeled effects in the integrated three-point function such as 
nonlinearities of the matter density, nonlocal bias parameters, and redshift-space distortions 
beyond the Kaiser factor, will tend to bias our measurement of cosmological parameters 
based on iC,{r). We will discuss caveats at the end of this section. 

Since the linear two-point and the tree-level three-point functions are proportional to (t| 
and (Tg, respectively, and is proportional to Ug, the scaling of the redshift-space correlation 
functions is 






(Tt 



O'S.M 


/bndry(^) 


(5.1) 

















baseline 

eff kernel 

tidal bias 

both 

b2 

0.41 ±0.41 

0.51 ±0.41 

0.48 ±0.41 

0.60 ±0.41 


Table 3: Best-fitting 62 and their uncertainties for BOSS DRIO CMASS sample for the 
extended models. The detailed description of the extended models is in appendix B. 


where “fid” denotes the quantities computed with the fiducial value of erg. Note that ^mc(^) 
is proportional to cr| because it is an integral of two linear power spectra (see eq. (3.9)). Since 
^/^o-„(r) dominates the signal, the parameter combinations bias and K = I + 2/3/3 -|- /3^/5 
are degenerate in the two-point function. That is, the amplitude of the two-point function 
measures only (^lUg)^ -|- |( 6 ict 8 )(/iT 8 ) -|- ^(/ug)^. This degeneracy can be lifted by including 
the quadrupole of the two-point function in redshift space. See refs. [37-40] for the latest 
measurements using the BOSS DRll sample. 

As for the three-point function, hgure 1 shows that the bf and 6^62 terms are comparable 
for 61 ~ 62 . This means that, at the three-point function level, the nonlinear bias appears in 
the leading order, so the amplitude of the three-point function measures a linear combination 
of bi and 62 - This provides a wonderful opportunity to determine 62 - The challenge is to 
break the degeneracy between 62 , bi, /, and ug. For this purpose, we combine our results 
with the two-point function in redshift space and the weak leasing measurements of BOSS 
galaxies. We take the constraints on bias{z = 0.57) = 1.29 ± 0.03 and f{z = 0.57)as{z = 
0.57) = 0.441 ± 0.043 from table 2 in ref. [37]. To further break the degeneracy between 61 , 
/, and as, we take the constraint on erg = 0.785 ± 0.044 from ref. [41, 42], where they jointly 
analyze the clustering and the galaxy-galaxy leasing using the BOSS DRll CMASS sample 
and the shape catalog from Canada France Hawaii Telescope Legacy Survey. 

We assume Gaussian priors on feicrg, fas, and ag with the known covariance between 
bias and fag- The cross-correlation coefficient between 6 icjg and fag is —0.59, as shown in 
figure 6 of ref. [37]. We then run the Markov Chain Monte Carlo with the Metropolis-Hastings 
algorithm to fit the model eq. (5.1) to the observed iC,{r')la\. We find 62 = 0.41 ± 0.41, and 
the results for the extended models are summarized in table 3. 

The value of 62 we find is lower than the mean of the mocks, 6 ™°^^ = 0.54 ± 0.35. The 
difference is mainly due to two reasons. First, the amplitude of the integrated three-point 
function of the observations is lower than that of the mocks by 10% {Ai = 0.89 ± 0.12). 
Second, the priors from the correlation function and lensing constraint bi to be close to 2.18, 
which is larger than that of the mocks, = 1.93. Thus, it requires a smaller 62 to fit the 

three-point function. The argument is similar for the extended models. Note, however, that 
the nonlinear bias of the data is still statistically consistent with the mocks. 

Let us conclude this section by listing three caveats regarding our cosmological inter¬ 
pretation of the measured integrated three-point function. 

1. The models we use, eq. (5.1), are based on tree-level perturbation theory, the lowest 
order redshift-space distortion treatment, as well as on the local bias parametrization. 
While this simple model describes the mocks well, as shown in section 3.3 and 3.4, we 
discuss in appendix B that using the effective F 2 and G 2 and the non-local tidal bias 
brings 62 closer to that of ref. [16]. We, however, find similar goodness of fit for various 
models, and thus we cannot distinguish between these models. 

2 . Covariances between the integrated three-point function, monopole and quadrupole 
two-point function, and weak lensing signals are ignored in our treatment. This can 
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and should be improved by performing a joint fit to all the observables. 

3. The cosmology is fixed throughout the analysis, except for / and a%. In principle, 
marginalizing over the cosmological parameters is necessary to obtain self-consistent 
results, although the normalized integrated three-point function is not sensitive to 
cosmological parameters such as flm as shown in hgure 6 of ref. [2] . 

These caveats need to be addressed in the future work. 

6 Conclusions 

In this paper, we have reported on the first measurement of the three-point function with 
the position-dependent correlation function from the SDSS-III BOSS DRIO CMASS sample. 
The correlation between the position-dependent correlation function measured within sub¬ 
volumes and the mean overdensities of those subvolumes is robustly detected at 7Aa. This 
correlation measures the integrated three-point function, which is the Fourier transform of 
the integrated bispectrum introduced in ref. [ 2 ], and is sensitive to the bispectrum in the 
squeezed configurations. 

Both the position-dependent correlation function and the mean overdensity are easier 
to measure than the three-point function. The computational expense for the two-point 
function is much cheaper than the three-point function estimator using the triplet-counting 
method. In addition, for a fixed size of the subvolume, the integrated three-point function 
depends only on one variable (i.e., separation), unlike the full three-point function which 
depends on three separations. This property allows for a useful compression of information 
in the three-point function in the squeezed configurations, and makes physical sense because 
the integrated three-point function measures how the small-scale two-point function, which 
depends only on the separation, responds to a long-wavelength fluctuation [ 2 ]. As there are 
only a small number of measurement bins, the covariance matrix of the integrated three- 
point function is easier to estimate than that of the full three-point function from a realistic 
number of mocks. We have demonstrated this advantage in the paper. 

Of course, since this technique measures the three-point function with one long-wavelength 
mode (mean overdensity in the subvolumes) and two relatively small-wavelength modes 
(position-dependent correlation function), it is not sensitive to the three-point function of 
other configurations, which were explored by ref. [16]. 

We have used the mock galaxy catalogs, which are constructed to match the two-point 
function of the SDSS-III BOSS DRIO CMASS sample in redshift space, to validate our 
method and theoretical model. We show that in both real and redshift space, the integrated 
three-point function of the mocks can be well described by the tree-level SPT model. However, 
the nonlinear bias which we obtain from the mocks is higher than that reported in ref. [16]. 
This is possibly due to the differences in the scales and conhgurations of the three-point 
function used for the analyses. As discussed in section 5, any unmodeled nonlinear effects 
in the redshift-space integrated three-point function of CMASS galaxies will tend to bias 62 , 
and will bias this parameter differently than the measurement of ref. [16]. 

Taking the mean of the mocks as the model, and treating the amplitudes of two- and 
three-point functions as free parameters, we find the best-fit amplitudes of iQ{r)/a\, i{r), 
and a\ of the CMASS sample. With respect to the mean of the mocks, the observations 
show a somewhat smaller {Ai = 0.89 ± 0.12) and larger a\, while the ensemble 

two-point function ^(r) matches the mocks. Given that the mocks are generated to match 
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specifically the two-point function of the BOSS DRIO CMASS sample within a certain range 
of separations, the level of agreement between the observations and mocks is satisfactory. 

Finally, by combining the integrated three-point function and the constraints from the 
anisotropic clustering {biag and fas in [37]) and from the weak leasing measurements (erg in 
[42]), we break the degeneracy between 6 i, 62 , f, and ug. We find 62 = 0.41 ± 0.41 for the 
BOSS DRIO CMASS sample. The caveat of this result is that our model, eq. (5.1), relies 
on a rather simple model in redshift space as well as on the local bias parametrization. We 
leave the extension of the model to improved bias and redshift-space distortion modeling 
(especially in light of the comparison with the results in ref. [16]) for future work. 

In summary, we have demonstrated that the integrated three-point function is a new 
observable which can be measured straightforwardly from galaxy surveys using basically the 
existing and routinely applied machinery to compute the two-point function, and has the 
potential to yield a useful constraint on the quadratic nonlinear bias parameter. Moreover, 
since the integrated three-point function is most sensitive to the bispectrum in the squeezed 
configurations, it is sensitive to primordial non-Gaussianity of the local type (parametrized 
by /nl)) thereby offering a probe of the physics of inflation. We plan to extend this work to 
search for the signature of primordial non-Gaussianity in the full BOSS galaxy sample. 
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A Testing the integrated three-point function estimator with Gaussian 

realizations and the local bias model 

We shall demonstrate that our integrated three-point function estimator is unbiased. To do 
this, we first generate the matter density field, Jm(r), by Gaussian realizations^ with Pi{k) at 

^Since Sm(r) follows the Gaussian statitistics, it is possible that Sm(r) < —1, which is unphysical. However, 
as we only compute the power spectrum without Poisson sampling the density field, this effect can be neglected. 
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Figure 8: The normalized integrated three-point functions of the mock halo density field with 
6 i = 3 and 62 = 1- The left and right panels are for Vl = 100 h~^ Mpc^ and 200 h~^ Mpc^, 
respectively. The data points show the mean of 300 Gaussian realizations, and the error bars 
are the variances of the mean (but note that the data points are highly correlated). The blue 
dashed and red solid lines are the theoretical models (iCfe) without and with [/bndry(^)]~^) 
respectively. 


z = 0 for the volume Vr of 1200 h ^ Mpc^ and a mesh size of 4 /i ^ Mpc. We then compute 
a mock “halo” density field using the local bias model via 


6hir) = bi6m{r) + y 




T 


rev; 


(A.l) 


where we set 61 = 3 and 62 = 1, and X^reF denotes a sum over grid cells in the entire volume. 
Note that X^rev;~ divide the entire volume W into Ng = 12^ = 1728 

subvolumes Vl of 100 h~^ Mpc^ and Ng = 6^ = 216 subvolumes Vl of 200 h~^ Mpc^. The 
two-point function in the subvolumes and the integrated three-point function are estimated 
by 


^hir,rL) = 


E 


x+r,xeVL 


4(x r)4(x) 


E 


x-|-r,xeVL 


1 


N, 

Kh{r) = Y, ih{'r,rL)6h{YL) , 

i=l 


(A.2) 


where S^i^L) is the mean halo overdensity in the subvolume centered at tl- Note that 
the denominator in the estimator of r^,) takes the boundary effect into account so 
r^,)) = without /bndry(^)- This means that the theoretical model of the integrated 
three-point function computed by eq. (2.9) has to be divided by /bndry(?’)- Since (5m (r) is 
Gaussian, the only contribution to the three-point function is from the nonlinear bias term, 
and so the estimated integrated three-point function is exactly given by 

Figure 8 shows the normalized integrated three-point functions of the mock halo density 
field with 61 = 3 and 62 = 1 from 300 Gaussian realizations. The measurements are in 
excellent agreement with gives us the confidence that our estimator is 

unbiased. 
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r-space 

hi 62 

baseline 
eff kernel 

tidal bias 
both 

1.971 ±0.076 0.58 ±0.31 

1.973 ±0.076 0.62 ±0.31 

1.971 ±0.076 0.64 ±0.31 

1.973 ±0.076 0.68 ±0.31 


z-space 

hi 62 

baseline 
eff kernel 

tidal bias 
both 

1.931 ±0.077 0.54 ±0.35 

1.933 ±0.077 0.65 ±0.35 

1.932 ±0.077 0.60 ±0.35 

1.933 ±0.077 0.71 ±0.35 


Table 4: Best-fitting values of hi and 62 and their uncertainties for mock catalogs, obtained 
using different models of the bispectrum in real space (left) and redshift space (right). 


B Effects of effective F 2 and G 2 kernels and non-local tidal bias 


In this appendix, we show how the inferred value of 62 changes when extending our baseline 
model for the bispectrum based on SPT at the tree level with local bias to the model used 
in the analysis of ref. [16]. 


Their model replaces F 2 and G 2 in eq. (2.21) with “effective” kernels, [34] and Gl® 
[35], which are calibrated to match the nonlinear matter bispectrum in of N-body simulations. 
Their model also adds a non-local galaxy bias caused by tidal fields [31-33] to Z 2 , i.e., 
Z 2 ^ Z 2 + ^bg 2 ^{ki ■ ^ 2 )^ — 5 , where 6^2 = —(4/7)(6i — 1). We use this model to compute 
the integrated three-point function, and find 62 of the mocks in real and redshift space by 
performing a joint ht with the two-point function as described in section 3.3 and 3.4. 


Table 4 summarizes the results. The “baseline model” refers to the model based on 
SPT and local bias. The “eff kernel” refers to the model with F|®, Gl®, and local bias. The 
“tidal bias” refers to the model with F 2 , G 2 , local bias, and tidal bias. Finally, “both” refers 
to the model with Ff®, G^®, local bias, and tidal bias. 

Both the effective kernels and the non-local tidal bias result in a larger nonlinear bias, 
which is in better agreement with ref. [16]. The changes of the best-fitting nonlinear bias, 
however, are still within the l—a uncertainties, and all the results are consistent with ref. [16]. 
We also calculate the goodness of the ht for all the models in both real and redshift space by 
comparing the mean of the mocks and the best-htting models, as well as the y^-distribution. 
We hnd that all models perform equally well; thus, in this paper we shall primarily use the 
simplest model, i.e. the SPT at the tree level with local bias for modeling the three-point 
function, but also report the results for the extended models. 


C Comparison of iC,(r)/a\ of BOSS DRIO CM ASS sample and PTHalos 
mock catalogs in different redshift bins 

The BOSS DRIO CMASS sample and the mocks have different sets of random samples 
with slightly different n{z), hence the properties of the observations and the mocks may not 
agree well in all redshift bins. Moreover, as mentioned in ref. [42], the CMASS sample is 
flux-limited, and thus the observed galaxies statistically have larger stellar masses at higher 
redshift (see figure 1 in ref. [42]). This may cause redshift evolution of the bias, and so the 
correlation functions. We shall measure i(^{r)/a\ as a function of redshift to test this. 

The measurements in the subvolumes are mostly the same as introduced in section 3.2, 
except that we now measure a{zj) as a function of redshift bin zj, and the average is done 
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Figure 9: iC,{r^la\ of 220 h ^ Mpc subvolumes in different redshift bins. The redshift bins 
increase from left to right, with the redshift cuts quoted in the beginning of section 3 and 4. 



<mock> 




30 40 50 60 70 80 

r [h- Mpc] 



Figure 10: Same as figure 9, but for 120 h ^ Mpc subvolumes. The redshift bins increase 
from top left to bottom right. 


in the individual redshift bin. Namely, 


a{zj) = 


E 

E 


i£zj ^ 9 ,* 




r,Zj 


i&Zj 


W, 




g{zj) = - 

Wr. - ^ 


i&Zj 


(C.l) 


This assures that 6{zj) = 0 for all redshift bins. 

Figure 9 and figure 10 show at different redshift bins for 220 and 120 h~^ Mpc 

subvolumes, respectively. We find no clear sign that iC,{r')la\ of the observations has different 
redshift evolution relative to the mocks. Thus, it is justified to study iQ{r)/a\ using one 
effective redshift for the BOSS DRIO CMASS sample. With the upcoming DR12 sample 
with a larger volume, the redshift evolution of iC,{r)ja‘j^ can be better studied. 


D Fisher matrix of the position-dependent power spectrum 

In this appendix we demonstrate that the position-dependent correlation function technique 
captures the information of the squeezed-limit bispectrum. To this end we utilize the lo- 
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cal type of primordial non-Gaussianity, which peaks in the squeezed limit, and show that 
the full bispectrum and the position-dependent correlation function would yield compara¬ 
ble constraints on /nl- As shown in section 2.3, the integrated three-point function is the 
Fourier transform of the integrated bispectrum. Here we shall compute the Fisher matrix 
of the integrated bispectrum instead of the integrated three-point function. This simplifies 
the computation as in this case we may assume the covariance matrix to be diagonal (in the 
Gaussian limit). 

As the main goal of this appendix is to show that the information of the squeezed-limit 
bispectrum is captured by the position-dependent correlation function technique, we shall 
use the simplest model for the galaxy bispectrum in the presence of local-type primordial 
non-Gaussianity and redshift-space distortion. Thus, the galaxy bispectrum is given by 

= Kiy [ 6 fi 3 spT(ki,k 2 ,k 3 ) -|- 6^62-862 (ki, k2, ks) -|- (ki, k2, ks)] , 

(D.l) 

where K^, = 1 -1- |/3 -|- is the Kaiser factor of the bispectrum and /? = f /hi [36], and 
/nl parametrizes the amplitude of the primordial non-Gaussianity. The exact forms of the 
bispectra are 


HspT(ki,k2,k3) = 2F2(ki,k2)P/(A:i,a)Pz(A:2,a) + 2 cyclic 
Hb2(ki,k2,k3) = Pi{ki,a)P 2 {k 2 ,a) +2 cyclic 
H/j^L(ki,k2,k3) = 2M{ki,a)M{k2,a)M{k3,a)[P<s>{ki)P,i.{k2)+ 2 cyclic] , (D.2) 


where a is the scale factor, M(k,a) = ^ k‘^T{k) and T{k) is the transfer function, and 

P^{k) is the power spectrum of the scalar potential. is produced by the local-type non- 

Gaussianity in the primordial scalar potential, <h(r) = </>(r)-|-/NL — (</>^(r))] where (/>(r) 

follows the Gaussian statistics [43]. Note that we ignore the effect of the scale-dependent 
bias due to the local-type primordial non-Gaussianity [44-46], and the more complete model 
can be found in [47, 48]. 

The Fisher matrix of the reduced bispectrum Q{ki, /c2, ks) = B{ki,k 2 , k 3 )/[P{ki)P{k 2 ) + 
2 cyclic] is given by [49] 


BQ,ai3 


E 

fcl ,^2 


dQjki, k2, ks) dQ{ki, k2, k^) _1_ 

dpa dpp AQ^{ki,k2,k3) 


(D.3) 


where {ki, ^2, k^) have to form a triangle, pa G [&i, 62 , /nl] are the parameters, and AQ‘^{ki, ^2, k^) 
is the variance of the reduced bispectrum. Similarly, the Fisher matrix of the normalized 
integrated bispectrum ib{k) = iB{k)/[P{k)(T‘j/\ is given by 


Fib,ai3 — 

L fc^fcmax 


dih{k) dih{k) 1 
dpa dpjs Aib‘^{k) ’ 


(D.4) 


where L is the subvolume size and Aib'^ik) is the variance of the normalized integrated bi¬ 
spectrum. Note that we assume the Gaussian limit in eqs. (D.3)-(D.4), hence only the diago¬ 
nal elements of Q and ib exist. Namely, for Q, it requires ki = k[ for {Q{ki, k 2 , k 3 )Q{k[, A:^, ^3)); 
similarly, for ib, it requires L = L' and k = k'. This approximation holds well at high redshift, 
but would break down at low redshift due to the nonlinear evolution. 
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Figure 11 : Two-dimensional joint 95% C.L. constraints on galaxy biases and primordial non- 
Gaussianity for BOSS full survey. The survey parameters are in the top right panel. The top- 
left, bottom-left, and bottom-right panels show the joint constraints on {bi, 62), {bi, /nl), and 
(^ 2 ) /nl) marginalized over /nl, ^ 2 , and bi , respectively. The blue dashed, green solid, and red 
dot-dashed line are for reduced bispectrum with kmm = 2TTjVr , reduced bispectrum with 
kmin = 2 vr/l/^' where Vl is the largest subvolume, and normalized integrated bispectrum for 
six sizes of subvolumes (100 to 600 h~^ Mpc with an increment of 100 h~^ Mpc), respectively. 


To compute the variances, we assume that the dominant component is from the bi¬ 
spectrum in the numerator (instead of the normalization in the denominator). We thus 
have 


AQ‘^{ki,k2,k3) 

Aib‘^{k) 


7rsi23 [Pzjkl) + Fshot] [Pz{k2) + Fshot] [^ 2 (^ 3 ) -I- Pshot] 
kik2k3 [Pz{ki)Pzik 2 y+ 2 cyclic ]2 

Vl kL + ^shot/VL][P.(fc) + .fshot] 


(D.5) 


where S 123 = 6 , 2, 1 for equilateral, isosceles, and general triangles, respectively, Pz{k) = 
Kpb\Pi{k) is the redshift-space power spectrum, Pshot is the shot noise of the power spectrum, 
NkL is the number of independent Fourier modes in the subvolume, and (y\ v 

Figure 11 shows the two-dimensional joint 95% C.L. constraints on galaxy biases and 
primordial non-Gaussianity for the full BOSS survey. One finds that as long as fcmin is 
set to be the fundamental frequency of the largest subvolume, the integrated bispectrum 
technique gives the similar constraint on /nl compared to the full bispectrum analysis.® 
This means that the integrated bispectrum is sensitive to the bispectrum in the squeezed- 

®Note that the numerical results are sensitive to the choices of fcmin and fcmax because we count the Fourier 
modes in this range. For different lines, though femax is set to be the same, in practice we stop counting Fourier 
modes if fc > fcmax- Therefore, they have different “true” fcmax, and the contour area would be affected. This 
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limit configurations. On the other hand, the top-left panel shows the strong degeneracy 
between hi and 62 for the integrated bispectrum. As discussed in section 2 . 2 , in the squeezed 
limit Bspt and Bb^ has the same scale-dependency if the power spectrum is pure power-law. 
Therefore, if only the squeezed-limit bispectrum is used, such as the integrated bispectrum 
technique, then it is challenging to break the degeneracy. 

We also find that while both techniques give similar constraints on /nL; the number 
of counted Fourier modes differ dramatically. Specifically, the full bispectrum using femin = 
271 jVr and 27r/y^' counts 7113 and 6730 triangles, respectively, while the integrated 
bispectrum counts only 54 Fourier modes. This means that many of the triangles contain 
information other than the squeezed configurations® and thus the integrated bispectrum 
can be regarded as an efficient approach to extract the information of the squeezed-limit 
bispectrum, such as the local-type primordial non-Gaussianity. 
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